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Help ed by the careful analysis of their experimental data, IWorthington fc Cold (|1897 . 

described roughly the mechanism underlying the formation of high-speed jets 
ejected after the impact of an axisymmetric solid on a liquid-air interface. They made 
the fundamental observation that the intensity of these sharp jets was intimately related 
to the formation of an axisymmetric air cavity in the wake of the impactor. In this work 
we combine detailed boundary-integral simulations with analytical modeling to describe 
the formation and break-up of such Worthington jets in two common physical systems: 
the impact of a circular disc on a liquid surface and the release of air bubbles from an 
underwater nozzle. We first show that the jet base dynamics can be predicted for both 
systems using our earlier model in Gekle, Gordillo, van der Meer and Lohse. Phys. Rev. 
Lett. 102 (2009). Nevertheless, our main point here is to present a model which allows us 
to accurately predict the shape of the entire jet. In our model, the flow structure inside 
the jet is divided into three different regions: The axial acceleration region, where the 
radial momentum of the incoming liquid is converted into axial momentum, the ballis- 
tic region, where fluid particles experience no further acceleration and move constantly 
with the velocity obtained at the end of the acceleration region and the jet tip region 
where the jet eventually breaks into droplets. Good agreement with numerics and some 
experimental data is found. Moreover, we find that, con trarily to the capillary breakup 
of liquid cylinders in vacuum studied by iRavleigh ( 18781) . the breakup of stretched liquid 



jets at high values of both Weber and Reynolds numbers is not triggered by the growth of 
perturbations coming from an external source of noise. Instead, the jet breaks up due to 
the capillary deceleration of the liquid at the tip which produces a corrugation to the jet 
shape. This perturbation, which is self-induced by the flow, will grow in time promoted 
by a capillary mechanism. Combining these three regions for the base, the jet, and the 
tip we are able to predict the exact shape evolution of Worthington jets ejected after the 
impact of a solid object - including the size of small droplets ejected from the tip due to 
a surface-tension driven instability - using as the single input parameters the minimum 
radius of the cavity and the flow field before the jet emerges. 



1. Introduction 

The impact of a solid object against a liquid interface is frequently accompanied by the 
ejection of a high speed jet emerging out of the liquid bulk into the air. Figure [lj which 
shows the effect of a horizontal disc that impacts on a pool of water, illustrates a liquid 
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Figure 1. Image of the high-speed jet ejected after the impact of 2cm disc with 1 m/s on a 

quiescent water surface. 



jet which flows ~ 20 times faster than the disc impact speed. The qualitative description 
of this common and striking phenomenon was firstly e lucidated at the beginning of the 
twentieth century bv lWorthington fc Cold (|l897l . llQODh . Through t he careful analysis of 
the photogr aphs taken after a solid sphere was dropped into water. IWorthington fc Cole 
( 18971 ll90Clh realized that these type of liquid threads emerge as a consequence of the 
hydro static collapse of the air-filled cav ity which is created at the wake of the impacting 
solid. IWorthington fc Cold (jl897l Il90fj ) also made the remarkable observation that the 
generation of such cavities was very much influenced by the surface proper ties of the 
spherical solid. One century after their original observations. iDuez et all (|2007l ) quantified 
the conditions that determine the existence of the air cavity in terms of the surface 
properties of the solid and the material properties of the liquid. 

High speed jets emerging out of a liquid interface are also frequently observed in 
many other situations. For instance, it is very usual to perceive that the liquid "jumps" 
out of the surface of sparkling drinks, a fact which is known to happen as a con - 



sequence of bubbles b u rsting at the liqu i d inte rface [Boulton-Stone fc Blakd (|1993I ); 
iDuchemin et al. ( 2002 ); Liger-Belair et al\ (120081)1 . Similarly , the impact of a drop o n 



a liquid int e rface or solid surface 
Reinl(|l993l ); lMorton et a7] (|2000h 



Qguz fc Prosperettil (fl99fjfi: IShin fc McMahonl (Il990h 



Deng et all (|2007l ) ; iBartolo et all (|2006l )]. is commonly 



accompanied by the ejection of liquid jets whose velocities can be substantially larger 
than that of the impa cting drop. Less fa miliar situations s u ch as t hose related to the 
focussing of capillary iMacIntvre dl968h : Thoroddsen et al. ( 20076h ] or Faraday waves 
[Hogrefe et all (|l998h ; IZeff et alT ( 200(h } also give rise to the same type of phenomenon. 
Nevertheless, in spite of the clear analogies, the main difference between the situa- 
tions enumerated above and the case of jet formation after cavity collapse is that, in 
the latter case, surface tension does n o t play any role in the jet ejection process [see 
Gekle. Gordillo. van der Meer fc Lohsd ( 2009al ) for details]. Indeed, the type of Wor- 



thington jets to be described here depend on a purely inertial mechanism, namely the 
radial energy focussing along the narrow cavity wall right before the cavity pinches-off. 
This fact makes our proce ss also somewhat d i fferent from situations in which jets are 
induced by pressure waves |Ohl fc Ikinkl (|2003h : lTjan fc Phillips! (|2007h : lAntkowiak et~al 
( 2007h : lBlake et all (|l993h ]. 
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Moreover, contr arily to what cou l d be e xpected from the analogy with ot her related 
physica l situations Longuet-Higgind (j 19831 ) : iLonguet-Higgins fc Qguzl(ll995l )]. lGekle etal 



( 2009a[ ) pointed out that jets formed after cavity collapse are not significantly influenced 
by the hyperbolic type of flow existing at the pinch-off location. Instead, the descrip- 
tion of this type of jets shares many similarities with the very violent jets of fl uidize d 
metal which are ejected after the explosion of lined cavities [e.g. Birkhoff et al. ( 19481 )]. 
with those formed whe n an axisymmetric bubble coll apses inside a stagnant liquid pool 
Manasseh et al. l (|l998l) ; [Boianos-J imenez et al. | ([20081 )] or possibly e ven with the granular 



j ets observed when a n obi ect impacts a fluidized granular material [Thoroddsen fe Shen 
(|200ll ); lLohse et al\ (|2004l )]. 

Most of the results presented here refer to the perpendicular impact of a circular disc 
with radius Rd and constant velocity Vd against a liquid surface. The fact that the solid 
is a disc instead of a sphere leads to the formation of an air cavity which is attached at the 
disc periphery, independent of the surface properties. Thus, this choice for the solid geom- 
etry avoids the additional difficulty of determining the position of the void attachment 
line on the solid surface. T he differences pointed out above set our system som ewhat 
apart from similar studies [Duclaux et al. ( 2007 ); Glasheen fc McMahonl (|l996t )]. The 
experimental realization of the setup to which the numerical simul a tions present ed are 
referred, is described by iBergmann et all (|2006l l2009h : iGekle et all (j2008l l2009af ). who 
show that boundary-integral simulations are in excellent agreement with experiments. In 
addition, potential flow numerical simulations to study of the type of Worthington jets 



water Manasseh et all (119981) :lLonguet-Higgins et al. (1991): Oguz & Prosoeretti 


dl993): 


iBurton et all (2005):lKeim et al. 


(l2006):lThoroddsen et all (2007a. 120081): Gordil 


et al. 


(I2007J); Burton & Taborek! ( 2008! 


):lGordillol(|2008l): Bolanos- Jimenez et al. (20081): 


Schmidt et al. 



)] are also reported in this paper. As in the case of Worthington jets ejected after 
solid body impact, similar boundary- integral simulations have b e en shown to be in excel- 
lent ag reement with experiments [see lOguz fc Prosperetti ( 19931 ) ; Bolanos- Jimenez et al 

ESS. 

This paper is organized as follows: In section[2]we present the three different numerical 
methods used. Section [3] presents the results from the simulations which are compared 
to the analytical model in section [4] Conclusions are drawn in section [5] 



2. Numerical methods 

In this paper we have used three types of boundary-integral simulations. The first two 
model, respectively, the normal impact of a disc on a free surface and the pinch-off of a 
bubble from an underwater nozzle. With the purpose of simulating the capillary breakup 
of the jets formed in the first two situations, the third type of simulation represents a jet 
issued from a constant-diameter nozzle with an imposed axial strain rate. The latter type 
of numerical simulations have the advantage of allowing us to directly impose the values 
of both the strain rate and the Weber number, which are the parameters controlling the 
breakup of the jet, as will become clear from the discussion below. 

2.1. Disc impact simulations 

The process of disc impact [see also IBergmann et al. (|2006ll2009l) : lGekle et al] (l2009al) ] 
is illustrated in figure [2] after impact a large cavity is created beneath the surface which 
subsequently collapses about halfway due to the hydrostatic pressure from the liquid 
bulk. From the closure location two high-speed jets are ejected up- and downwards. Here 
positions, velocities and time are made dimcnsionlcss using as characteristic quantities 
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the disc radius Rd, the impact velocity Vd, and To = Rd/Vd, respectively. (Variables 
in capital letters will be used to denote dimensional quantities whereas their lower case 
analogs will indicate the corresponding dimcnsionlcss variable). Moreover, it will be as- 
sumed that axisymmctry is preserved and, thus, a polar coordinate system (r, z) will be 
used. The origins of both the axial polar coordinate z and of time t are set at the cavity 
pinch-off height and at the pinch-off instant, respectively. 

Since gl obal and local R e ynold s numbers are large a nd the generation of vorticity is 
negligible [Bergmann et al. ( 20091 ); Gekle et al. (2009a)] we can make use of a flow po- 
tential to describe the liquid flow field. The numerical details, including the "surface 
surgery" needed to accurately capture th e transition f rom the cavity collapse process to 
the jet ejection, are given elsewhere [see Gekle et al. (2009a); Bergmann et al. ( 2009[ )]. 
These simulations have shown excellent agreement with experimenta l high- s peed record 
ings and part ic le image velocimetry me asurements Bergmann et al. ( 2006h : Gekle et al 
(120081 . 12009 of ): iBergmann et all (j2009h ] . The simulation stops when the downward jet 
hits the disc surface. 

Since the Reynolds number is large, the dimcnsionlcss parameters controlling the jet 
ejection process are the Froudc number, Fr = V^/ (Rd g), and the Weber number, 
Wc = pVpRo / & where g, p and a indicate the gravitational acceleration, the liquid 
density and the interfacial tension, respectively. Since Wc > O( 10 2 ) in all cases co n- 
sidered here, the jet ejection is not promoted by surface tension [Gekle et al. ( 2009a )]. 
Nonetheless, capillarity is essential to describe the jet breakup process, as will become 
clear from the discussion bel ow. Air effects, which play an essential role during the lat- 
est st a ges of cavity co l lapse IGordillo. Sevilla. Rodriguez-Rodriguez fe Martinez-Bazan 
(j2005t) ; lGordmol (|2008t ); lGekle. Peters. Gordillo. van der Meer fe Lohsd (|2009ij) ]. are not 
taken into explicit consideration here. Instead the cut-off radius at which the cavity ge- 
ometry is changed into the jet geometry is fixed manually verifying carefully that the 
exact value of this parameter does not influence our results. The only consequence of 
this simplification is that a tiny fraction of the jet - the jet tip - may not be accurately 
described neither by our numerical simulations nor by our theory as will be discussed in 
section 13.11 



2.2. Bubble pinch- off from an underwater nozzle 
In the second type of simulations a bubble grows and detaches when a c onstant gas flow 



rate is injec ted from an underwater nozzle into a quiescent pool of liquid. iManasseh et al 



( 1998| ) and iBolanos-Jimenez et al'. ( 2008 1 experimentally showed that this process also 
creates high speed jets. Indeed, as the bubble grows in size, the neck becomes more 
and more elongated and, eventually, surface tension triggers the pinch-off of the bubble, 
leading to the formation of two fast and small jets as illustrated in figure [3] Surface 
tension also leads to the pinch-off of a small droplet at the jet tip, which is precisely the 
instant when the simulation stops. 

Here, distances are made non-dimensional using the nozzle radius Rn as the charac- 
teristic length scale; moreover, the prescribed gas flow rate Q is used to derive the typical 
time scale T/v = (-kR? n )/Q. For the quasi-static injection conditions considered here, the 
relevant dimensionl ess parameter characte r izing t his physical situation is the Bo nd num- 



ber Bo = pR 2 N g/a [Longuet-Higgins et al. ( 1991 ); Bolanos- Jimenez et al. ( 20081) ]. which 



in the case presented here e quals 2.1. More deta ils of our simulation method are given in 
Qguz fe ProsperettH (Il993 ): Gekle et al. ( 2009c ). Note that the present numerical simula- 



tions 



method [Bolanos- Jimenez et all ( "20081 )]. are in excellent agreement with experiments 



uz fe Prosperetti j 19931)] as w ell as those reported using a very similar numerical 
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Figure 2. Numerical results obtained when a circular disc (blue line) impacts perpendicularly 
and at constant velocity on a flat liquid interface. Upon impact a cavity attached at the disc 
periphery is created in the liquid (a) which collapses under the influence of hydrostatic pressure 
(b). As a consequence of the cavity collapse, two jets with velocities much larger than that of 
the impact solid, are ejected upwards and downwards. The influence of increasing the impact 
Froude number from Fr = 5.1 - top row - to Fr = 92 - bottom row - is that the cavity becomes 
more slender. 
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Figure 3. (a) Time evolution of jets formed after the collapse of gas bubbles injected into a qui- 
escent liquid pool through a nozzle (red line), showing the ejection of the first drop, for Bo = 2.1. 
b) Closeup view of the jet region in (a). The colors correspond to different dimensionless times: 
t = (blue), t = 0.0014 (black) and t = 0.0027 (green) 
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Figure 4. Sketch defining the geometry of the numerical simulations used to describe the 
capillary breakup of a stretched liquid jet of density p injected into a gaseous atmosphere of 
density p g = 10 _3 p. The liquid velocity profile imposed at the boundary which delimits the 
nozzle on the left is uniform and decreases linearly with time. 

2.3. Simulations of a jet ejected at constant diameter 

As will be shown by our theoretical analysis below, the jet breakup process can be de- 
scribed in terms of two dimensionless parameters evaluated nearby the base of the jet, 
namely, the local Weber number and the dimensionless axial strain rate. These quanti- 
ties depend non-trivially on the input parameters of our physical simulations (disc speed, 
nozzle size etc.). In order to obtain a way of systematically varying both the local Weber 
number and strain rate we conducted a third type of simul ation by adapt i ng th e ax- 



isymmetric (two- fluid) boundary integral method described in I Gordillo et ali (|20071 ) to a 
situation that retains the essential ingredients to describe the capillary breakup process 
in the first two types of simulations. For this purpose, we have simulated the discharge 
of a liquid injected through a constant radius needle with a length of 20 times its radius 
into a gaseous atmosphere. The density ratio of the inner and outer fluids is 10 3 and a 
uniform velocity profile linearly decreasing with time is imposed on the boundary that 
delimits the computational domain on the left (see figure H|). Initially, the liquid interface 
is assumed to be a hemisphere attached at the nozzle tip. The uniform velocity with 
which the liquid is injected varies in time according to 

U N (t N ) = U N (0){l-at N ) (2.1) 

with the dimensionless strain rate a and the initial velocity Un(0) determined by the 
physical situation which one intends to imitate (jets formed cither after the disc impact 
or from the underwater nozzle). For these type of simulations positions, velocities and 
time will be made non dimensional using, as characteristic dimensional quantities, the 
injection needle radius R^, the initial velocity £/jv(0), and T N = R n /Un(0) respectively. 

In section 13.41 we demonstrate very good agreement between the results of these type 
of simulations and those related to the formation of jets after bubble pinch-off from 
an underwater nozzle. Unfortunately, the extremely large values of the Weber number 
reached at the tip of the liquid jets formed after the impact of a disc on a free surface (~ 
O(10 3 )) unavoidably lead to the development of numerical instabilities Tjan fc Phillips! 



< 20071 )]. This fact makes a direct comparison between the simulations of the axial strain 



system sketched in figure H] and those corresponding to the impacting disc impossible. 



3. Analysis of numerical results 

3.1. Effects of azimuthal asymmetries in the determination of the cut-off radius 

The value of r m i n (the minimum radius of the cavity before the jet emerges) would be zero 
under the ideal conditions of our simulations, which do not take into account gas effects 
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liquid viscosity [Burton et all 



iThoroddsen et all ()2007d ): iBolanos- Jimenez et 



(200S)1 or small azim u thal a symmetries that may be present in the flow 



Keim et 



2006); ISchmidt et al. (2009)]. This would imply that the initial jet velocity would be 



infinity. However, all the effects enumerated above are known to strongly influence the 
spatial region surrounding the cavity neck during the very last stages of bubble pinch- 
off and, therefore, are essential to determine the real value of r min (jGordilld (|2008f ); 



Gekle et al 



(1200961) 1. 

Note first that, the larger r m i„ is, the smaller will be the maximum liquid velocity 
at the tip of the jet. Here we will provide experimental evidence showing that non- 
axisymmetric perturbations are of crucial importance to fix r m in and, consequently, the 
maximum velocity reached by the jet. This is due to the fact that asymmetries influence 
the radial flow focussing effect on the central axis even before the actual cavity closure. 
The development of azimuthal instabilities leads to a decrease of the liquid acceleration 
towards the axis before pinch-off and thus reduces the speed of the ejected jet. This is 
clearly observed in figures [5] and [6l which show the cavity formation and jet ejection 
processes when either a brass disc (smooth surface) or a golf ball (structured surface) 
impact perpendicularly on a quiescent pool of water. Despite the fact that both the 
velocity and the diameter of the ball are larger than those of the disc, the maximum jet 
velocity is larger for the disc case. Indeed, while the shape of the cavity in figure [5] is 
smooth, the cavity interface in figure [6] clearly exhibits asymmetric modulations already 
right after the impact (which - in addition to the rough surface structure - may in part 
also be due to a rotation of the ball). Note that the overall shape of the cavity is very 
similar in both cases. Consequently, since the self-acceleration of the liquid towards the 
axis is lost when the amplitude of azimuthal disturbances is similar to the radius of 
the cavity, the maximum velocity reached during the collapse process decreases when 
the cavity interface is not smooth. Note that figures [5] and [6] are representative of an 
exhaustive set of experiments. The analysis of the whole experimental data has shown 
that the rough surface systematically produces lower jet speeds. 

The initial amplitude or the precise instant at which su c h azimuthal instabiliti e s ma; 



ine initial amplitude or tne precise instant at wmcn su c n azimutnai ins tabilitie s may 
develop is not easy to predict. For instance, iKeim et al. (|2006t ); ISchmidt et HI (|2009h 



pointed out that tiny geometrical asymmetries in the initial setup might break the cylin- 
drical symmetry of the cavity at the pinch-off location. Moreover, even if the cavity is 
perfectly axisymmetric, the strong shear between the gas and the li quid will induce in- 
stabili t ies that tend to break t he cylindrical symmetry of the cavity [Leppinen fc Lister 
(|2003t ); [BCTgmann et all (j2006l )]. 

Therefore, the precise determination of r m ;„ is a very complex and difficult subject 
which in addition will heavily depend on the system under study and must therefore 
remain outside the scope of this contribution. We have instead decided to vary r m i n within 
reasonable bounds and to analyze carefully the effect on the subsequent time evolution 
of the jet. It can be clearly appreciated in figure [7] that differences in the simulations can 
be observed in both the jet base and tip region right after pinch-off occurs. However, as 
soon as the jet radius at its base becomes of the order of the maximum value of r m i n 
explored, differences in the jet base region disappear and only remain appreciable in the 
jet tip region. Physically, this means that gas effects and small asymmetries will only be 
felt at the highest part of the jet, which represents only a very small fraction of both 
the total volume and of the total kinetic energy of the jet. Note also that, in spite of 
the jet tip being the spatial region where the highest velocities are reached, it is also the 
least reproducible one from an experimental point of view since it strongly depends on 
the precise details of pinch-off. Thus, regarding experimental reproducibility, our study 
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a:fl ms b:l>2 ms c. 64 ms 




FOnw fc 1.72 ms i: 4.09 ms j: 9.67 ms 



Figure 5. Pictures (a)-(f) show the smooth cavity formed after the normal impact of a brass 
disc against a water interface. The disc dimensions are 22 mm in diameter and 4.7 mm in height. 
The disc falls by gravity and the impact velocity is Vi mpa ct = 1-85 m/s. Note that, while the time 
between impact and cavity closure is roughly 70 ms, the upwards jet reaches the free surface 
in less than 4 ms, indicating that the jet velocity is much larger than the impactor's velocity. 
Indeed, the initial velocity of the tip of the jet, measured from detailed images of the type (g)-(j), 
is larger than - since drops might not be in a plane perpendicular to the free surface - 22.71 m/s 
and thus larger than 12.28 times the disc velocity. The huge velocities reached by the liquid jet 
can also be visually appreciated by comparison with the velocity of the drops formed in the 
corona splash which hardly change their position between images (g) and (j). Let us also remark 
that, initially, the jet is not axisymmetric ((h) and (i)). Nevertheless, after a few milliseconds, 
picture (j) shows that the jet becomes approximately axisymmetric. 
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Figure 6. Pictures (a)-(f) show the cavity formation caused by a golf ball with a diameter 
of 42.75 mm impacting with a velocity of 2.03 m/s. Compared to figure [5] the surface shape 
is visibly distorted (c) due to the rough surface structure of the ball. Nevertheless, it can be 
inferred from a detailed image analysis that the jet velocity is again much larger than the ball's 
velocity. However, in spite of both the impact velocity and the ball diameter being larger than 
those of the disc, the maximum velocity of the jet is only VSmpact — 20 m/s and thus smaller 
than for the impacting disc. 



will be valid to accurately describe the most robust part of the jet. In the case of the 
impacting disc we will set r min = 0.01 and in the case of the gas injection needle, the 



minimum radius will be fixed to tv 



0.05. 



Finally, note that our axisymmetric approach has been proven to be in excellent agree- 
ment with experiments whenever either the radius of the collapsing cavity or the radius 
of the emerging jet, are larger than the cut-off radius r m i n for which any of the effects 
enumerated above - gas , azimuthal pertur b ations - become relevant [see , for instance , 



(|2009cJ ); lBolanos-Jimenez et all (120091 )] 



Bergmann et al\ (120061); iGekle et all (20081 ); iBolanos- Jimenez et al\ (|2008l ); iGekle et al 



3.2. Jet ejection process for the disc impact 

The different stages of the jet formation process have been illustrated in figure [21 After 
the solid body impacts against the free surface, an air cavity is generated (a). As a 
consequence of the favorable pressure gradient existing from the bulk of the liquid to the 
cavity interface, the liquid is accelerated inwards (b). These radially inward velocities 
focus the liquid towards the axis of symmetry, leading to the formation of two fast and 
sharp fluid jets shooting up- and downwards, as depicted in figured] (c). Here we will 
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5.1 at two different instants of time, 



3 (b), for four different values of the cut-off radius, r m i n — 0.005 
= 0.02 (black), and r m i n = 0.05 (green). Bottom row: Jet shapes 



Figure 7. Top row: Jet shapes for the disc impact at Fr 
t = 10" 4 (a) and t = 3.2 x 10 
(blue), r min = 0.01 (red), r min 

for the underwater nozzle with different cut-off radii (colors as in top row) at t — 0.0003 (c) and 
t = 0.004 (d), respectively (here the simulations are extended beyond the ejection of the first 
droplet). It is evident in both cases that the influence of varying the cut-off is significant only 
in the very first instants after pinch-off and at the very tip of the jet. 



mainly focus on the detailed description of the upwards jet and demonstrate that the 
downward jet can be treated in the same way. 

From figure [H observe that larger Froude numbers create more slender cavities and 
also increase the non-dimensional depth at which the cavity pinches-off. Furthermore, it 
can be appreciated that the jets are extremely thin and that the time needed for the tip 
of the jet to reach the free surface is only a small fraction of the pinch-off time. This 
latter observation means that the jets possess a much faster velocity than the velocity 
of the impacting solid, a conclusion which was also extracted from the analysis of the 
experiments in figures 0(51 Motivated by this striking fact, one of the main objectives in 
this paper will be to address the following question: what is the relationship between the 
impact velocity Vd - or, in dimcnsionlcss terms, between the Froude number - and the 
liquid velocity within the jet? 

With this purpose in mind, it will prove conveni ent to define first the length scale that 
characterizes the jet width. In lGekle et al. ( 2009a ') we showed that the time evolution of 
the jet is a local phenomenon, independent of the stagnation-point type of flow generated 
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Figure 8. Sketch showing the different lengths used to define the jet base and the regions of 
the jet. The jet base (rt, Zb) is located where the interface possesses a local minimum. The outer 
region covers the bulk of the fluid with r > r;, and z < z\>. The jet region is subdivided into the 
acceleration, the ballistic and the tip region. Note that, in the following, it and v will be used 
to denote axial and radial velocities, respectively 



after pinch-off at the location where the cavity collapses. Therefore, this c haracteristic 



length n eeds to be related to a local instead of a global quantity and, following lGekle et 



( 2009 ah . we choose the radial position at which the interface possesses a local minimum 
i.e., the radius rb{t) indicated in figure [H We shall in the following call this point the jet 
base and denote its vertical position by Zb(t). 

To clearly show the spatial region surrounding the jet base, some of the different jet 
shapes taken from the time evolutions of figure O arc translated vertically so that they 
share a common vertical origin, as depicted in figure O Note that both the jet base and 
the jet itself widen as the time from pinch-off increases. Interestingly enough, figure [10] 
shows that jet shapes exhibit some degree of self-similarity since they nearly collapse onto 
the same curve when distances are normalized using r^. This fact indicates that rt, is not 
an arbitrary choice, but a relevant local length that plays a key role in the dynamics of 
the jet. The same arguments hold for the downward jet as illustrated in figure [Til 

In order to model the full process of jet ejection and break-up we divide the liquid 
flow field into two different regions: the outer region, defined for r > r^, z < Zb and the 
jet region, extending from the jet base to the axis i.e, r < and z ^ Z&, as illustrated 
in figure [TUl The jet region is further divided into three different axial subregions: the 
acceleration region, the ballistic region and the tip region as illustrated in figure [8l 

Figure [P2l shows that Zb(t) ^> rb(t). These comparatively large values of Zb{t) with 
respect to r^it) are caused by the confinement of the jet by the cavity walls, which 
inhibits the widening of the base radius. Moreover, the small values of are responsible 
for the large axial velocities within the jet (and, thus, for the large values of dzb/dt) 
since, as it will become clear below, vertical velocities are inversely proportional to rt,. 

The importance of local processes around the jet base is even more clearly illustrated 
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Figure 9. Jet shapes translated vertically for different instants of time and different values of 

the impact Froude number. 




Figure f 0. Shapes of the jets depicted in figure[9]when distances are normalized using r;, overlay 
reasonably well indicating that rt is a good choice for the characteristic local length scale. The 
definition of the outer and jet regions is also indicated. 
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Figure 1 f . Jet shapes for the downward jet at Fr=92 taken at the same times as in figures [9] 
and [10] show a similarly good rescaling with the corresponding base radius r^. 
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t t 

Figure 12. Time evolution of radial and axial positions of the jet base, r& and Zb respectively. 
The upward jet is shown in black and the downward in red (for the downward jet — Zb is shown 
for convenience). The behavior of both jets is very similar. 

in figure [13] where both the axial and radial velocities evaluated at the jet air/liquid 
interface, u and v respectively, are represented for different instants of time. In this 
figure one can observe that while the axial velocities are of similar magnitude as the radial 
velocities at r = r b , they monotonically increase to much higher values as the jet radius 
diminishes. Contrarily, the modulus of the (negative) radial velocities decays from ~ 
O(10) at r = Tb to zero at r ~ 0.5 and, therefore, the radial inflow experiences a strong 
deceleration in the small distance ~ 0.5r;,. Since the liquid is at atmospheric pressure at 
the free surface of the jet, the strong radial deceleration provokes an overpressure below 
the jet base. Accordingly, a strong favorable vertical pressure gradient is created and, 
therefore, the liquid experiences a large upwards acceleration in the vertical direction, 
creating the high speed jet ejected into the atmosphere. 

In the following, we shall define ro = 0.5rj, as the radial position on the jet interface 
at which radial velocities become negligible -v « for jet radii smaller than ro- and 
the corresponding vertical position and velocity, will be denoted in what follows Zq and 
Uo, respectively. Moreover, we will also define at this point a local Weber number as 
Weo = p Uq Ro ja = We Uq r whose time evolution is depicted in figure Q3J The large 
values indicate that surface tension effects can be neglected in the description of the jet 
ejection process. 

Thus, since the jet interface can be considered to be at constant atmospheric pres- 
sure and surface tension effects are negligible near the jet base, the only source of axial 
acceleration is the axial pressure gradient caused by the radial deceleration of the flow. 
Remarkably, this radial deceleration takes place in a very localized region nearby the jet 
base. (For radial positions on the jet smaller than ro already v ~ as shown in figure 
[131 ) Therefore, the source of axial acceleration (radial deceleration) is no longer active 
high up into the jet, but only near the jet base. This key observation is used to define two 
of the three different regions within the jet: the axial acceleration region for ro < r < r& 
and Zb < z < zq and the ballistic region for r < ro, z zq. The term used to name the 
latter region is based on the fact that, since v ~ for r < ro and the pressure at the jet 
interface is atmospheric, the momentum equation projected in the axial direction yields 

D u 

-j^j- = for z > zg with !i^«(r), (3-1) 
and D / Dt indicating the material derivative. Equation (|3.ip implies that fluid particles 





Figure 14. Time evolution ol the local Weber number at the beginning of the ballistic region 
for two different values of the impact Froude number. The large values demonstrate that surface 
tension is not relevant during the jet ejection process. To facilitate the comparison between 
the different Froude cases, times have been normalized by t ma x, which is the time when the 
downward jet hits the disc and the simulation stops. 



are no longer accelerated upwards and conserve the vertical velocities they possess at 
z = zq, which is the axial boundary between the axial acceleration region and the ballistic 
region. In equation (|3.ip . u ^ u(r) since the radial velocity gradients of axial velocities 
are negligible in the ballistic region (not shown). 

As a next step, we would like to scale the radial velocity field in the vicinity of the 
jet base which is, as discussed above, the source of momentum driving the jet ejection 
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Figure 15. The sink strength distribution q c (z) for two values of the Froude number. 




Figure 16. Spatial and temporal evolutions of the radial and axial velocities depicted in figure 
1131 normalized with qb(t)/rt(t). Observe that u/(qb/ri,) and v/(qt,/ri,) nearly collapse onto the 
same curves for each of the two values of the impact Froude number considered, Fr = 5.1 (a) 
and Fr = 92 (b). 



process. These radially inward velocities are originally created by the difference between 
the hydrostatic pressure in the bulk of the liquid and the gas pressure inside the cavity. 
After pinch-off however, the radial velocity field feeding the liquid jet is not appreciably 
modified by gravity during the time evolution of the jet since the local Froude number at 
the beginning of the ballistic region is Fro = U§ / (gZ sm f acc ) 3> 1 with Z sur face the axial 
distance between the beginning of the ballistic region and the height of the free surface 
far from the impact region (see figure [2j). 

Therefore, the radial velocities which give rise to the jet emergence can be charac- 
terized by the sink strength distribution at t = right before pinch-off occurs: q c (z) = 
— r c (z)r c (z), where r c an d r r indicate the rad ius of the cavity and its associated radial 
velocity, respectively (see Gekle et al. ( 2009 oh ). The values of q c (z) are shown in figure 
1151 for several values of the impact Froude number. 

In order to demonstrate the intimate relation of the jet ejection process with the 
velocity field right before pinch-off, we normalize the velocities v and u at the jet surface 
(as shown in figure IT5|) using, as the characteristic scale for velocities, qb(t)/rb{t), where 
qb(t) = q c {z = Zb(t)) is the sink strength at the height of the jet base. The remarkable 
result, depicted in figure HU is that both rescaled velocities nearly collapse onto the same 
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Figure 17. (a) The time evolution of Bt demonstrates that Bt is roughly constant in time, but 
does depend somewhat on the Froude number, (b) Taking the average of Bt over time (with 
error bars indicating the min/max) for different Froude numbers yields a function _B(Fr) which 
varies only between 2.5 and 3.5 in the range 3 ^ Fr ^ 92. As indicated in figure [141 i max is the 
time when the downward jet hits the disc and the simulation stops. 



master curves for a given Froude number and thus are almost constant in time for a 
fixed value of the rescaled position r/ri, < 1. This implies that, for a fixed value of 
axial velocities are inversely proportional to ri, i.e., the smaller the jet base radius - or, 
cquivalcntly, the more confined is the jet by the cavity walls -, the larger will be the axial 
liquid velocities within the jet. 

Of critical importance for our forthcoming discussion is the rescaled axial velocity 
evaluated at the boundary of the ballistic region, B t = u (t)/(qt(t)/rb{t)), whose time 
evolution is depicted in figure [17] (a). In accordance with the collapse of the rescaled 
velocities on a single master curve depicted in figure 1161 B t hardly changes with time 
and, thus, we can define the function B(Fr) = Mo/(9b/ r f>) which depends also very weakly 
on the Froude number, as depicted in figure [T71 (b). 

The result in figure [17] possesses the additional remarkable implication that axial ve- 
locities within the jet are larger than the radial velocities existing at the cavity inter- 
face before pinch-off occurs. This can be seen directly by recalling that \qb/n\ = \fb\, 
such that B is the ratio between the axial velocity uq with which fluid is ejected into 
the jet and the radial inward velocity at the jet base. Then, during the initial in- 
stants of jet formation, ~ r m i n , with r m j„ the minimum radius of the cavity be- 
fore jet emerges. Therefore, since the maximum radial velocity before pinch-off occurs 
is l^mml = \Qc(z = 0)/r m i„|, the maximum axial velocity within the jet is given by 
max(ito) = B(Fr)q c (z = 0)/r m i n ~ 3r m i„. This means that, essentially, the velocity 
with which the jet is ejected is roughly three times larger than the maximum radial 
velocity attained before pinch-off! 

In addition, provided that Weo » 1, fluid particles conserve their axial velocity within 
the ballistic region [see equation (j3.ip ] and, consequently, the tip of the jet transports away 
from the pinch-off location very valuable information about the largest velocities reached 
during the cavity collapse process. The knowledge of the function B could thus allow an 
experimentalist to estimate the maximal pinch-off velocity simply from measurements of 
the jet tip velocity. 
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Figure 18. The time evolution of the shapes of the jets ejected after bubble pinch-off from an 
underwater nozzle, show good overlap when distances are normalized using tv 



3.3. Jet ejection after bubble pinch- off from an underwater nozzle 

This section is devoted to the s tudy of Worthington j ets which are ejected after th e 
bubble collapse into a liquid pool Manasseh et al. ( 1998 ); Bolanos- Jimenez et al. ( 2008f )]. 
As depicted in figure [31 these jets are quite similar to the ones formed after the impact 
of a solid body against a free surface and, thus, we expect that the conclusions of section 
13.21 can be also used for their description. 

Figure [181 shows that, similarly to section [3~2l the different shapes nearly collapse onto 
the same master curve when distances are normalized using r&. This fact corroborates 
that ri, is also the correct length scale to characterize this type of jets. However, differently 
to the case of Worthington jets ejected after the impact of a solid body against a free 
surface, in which Weo > 10 3 , the local Weber number evaluated at the beginning of the 
ballistic region is ~ O(10) in this case (see figure fT9k). As a consequence of this, the total 
length at breakup of these jets is ~ 0(r\,) (see figure ITS"]), i.e, much shorter than the length 
of the Worthington jets in section 13.21 Moreover, such comparatively low values of the 
local Weber number indicate that surface tension has an effect in the description of the 
jet ballistic region. This is clearly appreciated in figures IT51 and l2"ul where the collapse onto 
each other of the normalized time evolutions of the axial and radial velocity components 
evaluated at the free surface {u and v) is also a bit deteriorated when compared with the 
case depicted in figure [TBI Nevertheless, the two main prerequisites for the model to be 
presented in the following sections are also satisfied in this case: firstly, the acceleration 
and ballistic regions are clearly differentiated in figure [20] and, secondly, the normalization 
of the interfacial velocities with q\, / r\, leads to a reasonable collapse onto a single master 
curve (see figure [20j) . 

3.4. Jet breakup 

The growth of capillary perturbations in a cylindri cal liquid jet, firstly quantified by 



ine growtn oi capillary perturbations m a cylindri cal liquid jet, nrstly quantmed Dy 
iRavleigbl (|l878h [see also lEggers fc Villerma"u"xl (|2008l )] , is based on the assumption that 
fluid particles conserve, to first order, their longitudinal velocity U. Rayleigh's analysis 
shows that, moving in a frame of reference with the jet velocity U (which in his case is 
constant U ^ U(z, t)), and no matter how large the Weber number is, the jet breaks due 
to the growth of capillary perturbations of wavelengths larger than the jet perimeter. 
The characteristic time needed for such perturbations to disrupt the jet into drops is the 
capillary time, ~ (pi? 3 /^) 1 / 2 , with R the jet radius. Therefore the jet breakup length, Lb, 
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Figure 19. (a) Time evolution of the local Weber number for the jet depicted in figure [3] (b) 
Time evolution of the normalized strain rate so,i oca i at the beginning of the ballistic region for 
the same case as in (a). 




Figure 20. (a) Axial and radial velocities evaluated at the jet interface for the case depicted 
in figure [3] In analogy with figure 1161 both the acceleration and ballistic regions are clearly 
identified: the modulus of the radial velocities decreases from r = to become negligibly small 
for r/rb < 0.5. (b) The same as in (a), but with velocities scaled with qb/rb. Due to the fact 
that the Weber number is considerably smaller in this case than for the impacting disc, the jet 
tip region can be appreciated in this figure as the multivalued part of the curves u and v for 
r/r b ~ 0.4. 



is such that Lb/ Ret (pU 2 R/a) 1 ^ 2 if aer odynamic effects are absent [Sterling fc Sleicher 
( 1975h ; Gordillo fc Perez-Saborid ( 20051 )]. Notice that the study of jet breakup in our 
case is somewhat related to that considered by Rayleigh since the fluid particles conserve 
their velocities, in a first approach, a long the ballisti c region of the jet. 

Similarly to the case considered bv lRavleigh (1878), the study of the capillary breakup 
of stretched jets will be divided in two: the calculation of the basic flow, which is free of 
capillary effects and the analysis of capillary waves propagating and growing in amplitude 
at the jet tip region. Viscous effects will be neglected in the analysis. 



3.4.1. Unperturbed flow 

If Wc 3> 1, the time evolution of both the jet radius and the liquid velocities in the 
ballistic region can be calculated neglecting surface tension forces an d making use of the 
slenderness of the jet. However, differently to the case considered by iRavleigh (1878), in 
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which the jet radius rj = 1 and u = constant, here u and rj are functions of z and t. In 
effect, if the fluid is assumed to follow a purely vertical motion inside the ballistic region, 
the couple of equations that determine u and rj are the momentum equation (|3.ip , which 
can be also written as 

Du du du ,„ „. 

7» =0 ^* + "&= ' (3 ' 2) 

and the unidirectional version of the continuity equation, namely, 

% + d -pl = 0^^1 = -p, (3.3, 
at Dt oz 

where D / Dt = d/dt + ud/dz indicates again the material derivative. From equations 
(|3.2j) - (|3.3j) . u, rj and Zj - the height at which the jet radius is rj - are completely deter- 
mined if the relevant quantities at the beginning of the ballistic region (tq, zq, and the 
velocity u ) are known functions of time. Indeed, equation (|3.2p expresses that fluid par- 
ticles conserve the vertical velocity they possess at the beginning of the ballistic region. 
Consequently, a particle ejected from the acceleration into the ballistic region at time 
t < t will, at time t, have attained a height 

Zj (t) =*b(r) + (t-r) «o(r). (3.4) 

To obtain the corresponding jet radius rj, equation (|3 ,3|) can be readily integrated to 
give 

r 2 (z -z-t)- r 2 (r) u (r)-dz /dT 

* 1 °' } ~ o[ ' u {t) - dz /dT - du (T)/dT(t - t) ■ [ ' 

Introducing the definition of the strain rate at the beginning of the ballistic region 

So r = -z-{z = zq) = rl r-rr (3.6) 

oz u (t) - z [t) 

allows us to rewrite equation (|3.5p in a more compact form as 

Note that Rayleigh's original analysis, u = constant and rj = 1 (cylindrical jet) may be 
recovered from equations (|3.2j) - (|3.3p by setting sq = and Uq, zq and ro constants in 
time. However, in our case, uq < and uo = Bq^/r^ > Zb — Zq and, therefore, by virtue 
of equation (|3.6p . these conditions imply s Q > 0; consequently, from equation (|3.7p . the 
jet is not cylindrical since it stretches downstream. 

Now, in order to obtain the complete jet shape at time t, we vary r between and 
t and use equations (|3.4|) and (|3-5[) to compute the corresponding vertical and radial 
coordinates of the jet. Note that, clearly, the particle ejected at r = will end up forming 
the tip of the jet. The comparison between the numerical results and those obtained 
from the integration of (|3.2p - (|3.3p . with the values of u (t), z d (t) and r Q (r) taken from 
the numerical simulations, is depicted in figure [^U The excellent agreement between 
numerics and the model validates the approach of considering that fluid particles conserve 
their axial velocities within the ballistic region. It should be pointed out, however, that 
equations (|3.2p - (|3.3p need to be corrected at the tip of the jet, where surface tension 
effects need to be retained. 

3.4.2. Growth of capillary disturbances 

The linear stability analysis for the type of velocity field given in section 13.4.11 was 



firstly accomplished bv lFrankel fc Weihsl (Il985h . who recovered Rayleigh's original result 




Figure 21. Comparison between the numerical jet shape and that obtained from equation (|3,5[1 
for the disc impact at Fr=5.1 (a) and the jet ejected from the underwater nozzle (b). The black 
line is the simulation and the red line is the analytical model. The input values of ro(t), zo(t) 
and uo(t) for the jet stretching model are taken from the simulations. Note that, since surface 
tension is not included in this model, the tip of the jet requires a separate treatment as described 
in section [3.4.21 



in the limit of s g = 0. It is our purpose here to extend the analysis on the breakup of 
stretched jets of Frankel fc Weihd (1985) to account for non linear effects and also for 
the influence of the tip. It is worthy to mention that, in our numerical approach, the 
wavelength of fastest growth rate is naturally selected by the local flo w around the jet 
tip an d, therefore, a linear stability analysis of the type reported in Frankel fc Weihsl 
(| 19851 ). is avoided. 

As a first step, the dimensional counterparts of r (r) [i?o(j)] and u (t) [Uo(t)] are 
chosen as the characteristic scales for lengths and velocities, respectively. Consequently, 
dimensional analysis indicates that the evolution of capillary perturbations in the ballistic 
region for t > r will solely depend on the dimensionless parameters Weo and soj oca i = 
S {t) R /U = s Q r Q /u Q . 

The values of Weo and so^ oca i depend non-trivially on the dimensionless parameters 
controlling the two different physical situations analyzed here. Consequently, in order to 
study systematically the jet breakup process as a function of So,iocaZ and Weo we employ 
the third type of simulations of the axial strain type described in section [2] and illustrated 
in figureSJ The real jet breakup process can then be reproduced provided that the values 
of the Weber number and the strain rate at the nozzle exit coincide with those at the 
beginning of the ballistic region i.e, Weo = We at = p 1/^(0) Rn/ct and so, local = a. 

Note that, since the values of the Weber number based on Un and p g are always such 
that We g = pglljj Rn / <J <C 1, the gas dynamics can be neglected and the only two 
relevant dimensionless parameters characterizing the axial strain system of figure 2] are 
Wejv and a. 

The numerical results depicted in figure [22] show a slender liquid thread which breaks 
many diameters downstream the nozzle exit. Moreover, it can be observed that the effect 
of increasing the Weber number is to increase the breakup time and the breakup length. 
Figure |2"51 shows a comparison between the shapes of the jets formed after bubble collapse 
depicted in figure [18] and those obtained from the simulations of the type illustrated in 
figure [22l with We^ = Wco(r = 0) and a = so t i OC ai(T = 0). The excellent agreement 
between both type of numerical results corroborates the fact that tip breakup of Wor- 
thington jets can be reproduced by means of the simulations considered in this section 
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Figure 22. Time evolution of jets calculated for two different values of the pair of variables 
(Wejv, a) but the same value of the product WejvQ 2 . 



(z-Ztip)/rdrop 1 




Figure 23. Jet shapes of figure [3] (in continuous lines) at three different instants of time 
(0.6T(, rea fc , 0.8Ttreak and Tbreak, with Ttreak the breakup time) compared to those obtained 
from the type of simulations depicted in figure [22] (dashed lines) for the same values of T. The 
values of the parameters have been set to Wejv = 25 and a = 0.4 which are the initial values 
of the local Weber number and the dimensionless axial strain rate of the simulations depicted 
in figure [19] Note that distances are rescaled using the final radius of the drop Rdrop as the 
characteristic length scale and that zu P denotes the axial coordinate of the tip of the jet. 



if the values of the parameter We at and a coincide with the initial values of Weo and 

SQ.local • 

However, the numerical code used in this section is unstable for We^v > O(10 3 ). Conse- 
quently, the third type of simulations cannot reproduce, at first sight, the breakup of jets 
ejected by an impacting disc since, in this case, Weo ^ 10 3 as depicted in figurefTll Thus, 
is it nevertheless possible to describe the breakup process of jets with such high values of 
Weo using the numerical simulations of the axial strain type illustrated in figure 2]? The 
answer to this question is affirmative if we realize that, in a frame of reference moving 
at the tip velocity, the parametrical dependence on the velocity Uq(t = 0) disappears. 
Consequently, since both the local flow field and the jet radius still depends in this frame 
of reference on Sq(t = 0) (see equation (|3.6p ). dimensional analysis indicates that jet 
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Figure 24. Translated jet shapes corresponding to the conditions Wejv = 90, a — 0.33 
and We at = 200 and a = 0.223, depicted in figures [22] (a) and (b), respectively. Since 
Weg = WejvQ 2 = 10 in both cases, the time evolution near the tip region is identical in the 
normalized temporal variable TSo- 



breakup can be described in terms of the dimcnsionlcss variables T S a (or, analogously, 



t N a) and Wc^ = p Sq Rq/ct = We rg s§ = Wc s 



.local 



(or, analogously, Weaver). 



To check this, the results depicted in figure [221 which correspond to different values of 
We 7v and a but to the same value of We at a 2 , are represented in figure |2~H Remarkably, 
the different jet shapes superimpose onto each other for the same values of the dimen- 
sionlcss time TSo, what indicates that the breakup process depends solely on Wes (or, 
equivalently, on Wejyck 2 ) and on the dimcnsionlcss time T So (or, equivalently, on ijva) 
for sufficiently large values of Weo. Figure [231 illustrates that the volume of the nearly 
spherical drops formed at breakup, decreases for increasing values of Wes- In figure [231 
note also that the range of values of Weg investigated is realistic even for the impacting 
disc, as depicted in[26] Consequently, even though Weo in some situations such as the disc 
may be very high, the important parameter which is Wes = Wcgs 2 , local can be matched 
to the simulations in this section. The dimcnsionlcss breakup time and the dimension- 



2/7 



less size of the drops in figure [57] behave as rdrop <x We s and (TS )break oc We 
respecti vely. A d e tailed corresponding theory will be the subject of a forthcoming contri- 
bution jGordillol (|2009t) ]. We emphasize that figure [57] describes a universal relation for 
the breakup of Worthington jets at high Weber numbers which allows one to obtain the 
breakup time and volume of the first ejected droplet knowing merely the value of We5 
defined at the beginning of jet formation. 

Finally, note that, in order for We at < O(10 3 ), the computations have been performed 
choosing 1 <C We at <C We and, correspondingly, a > sojocai- The condition We at 3> 
1 is essential since, if WeAr was not sufficiently large, the jet breakup process of real 
Worthington jets would depend on the liquid velocity Un (0) and, thus, on WeAr and a 
separately. 



4. Modeling the jet ejection and breakup processes 

Here we aim to develop a model to explain, in simple terms, the jet ejection and 
breakup processes. Our model will be based on the main conclusions of the previous 
section which arc: (i) both r& and Zb are local quantities which, therefore, do not depend 
on the large scales of the flow, (ii) the velocity field within the jet can be characterized 
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Figure 25. Tip region at the instant of breakup for different values of Weg. Observe that the 
drops generated are nearly spherical and that their volume decreases for increasing values of 
We s . 
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Figure 26. Time evolution of the normalized Weber number Wes = Weo Sq iocai for the im- 
pacting disc. As indicated in figure ITU tmax is the time when the downward jet hits the disc and 
the simulation stops. 



(TS )break 




Figure 27. Numerical results of the dimensionless breakup time, (TSo)break (a), and the 
dimensionless equivalent radius, r eq (b), with A/Sirr 'iq Rq the volume of the drops formed. 



solely in terms of the sink strength intensity at pinch-off, q c (z) and (iii) the flow field 
within the jet can be divided in three parts: the acceleration region, the ballistic region 
and the jet tip region. 

This section is structured as follows: in subsection 14. II n, (t) and Zb(t) are calculated in 
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Figure 28. The sink distribution q c {z) at the moment of pinch-off (green) for the impacting disc 
at Fr=5.1 is the essential ingredient to our jet formation model. The sink distributions at later 
times (red, blue and black curves) are almost unchanged with respect to the sink distribution 
at pinch-off, confirming our model assumption that q c (z) is valid even for jet formation when 
two additional effects are accounted for: the accumulation of sinks around the base and the hole 
between the up- and the downward jet. 



terms of only q c (z) using the theory developed in lGekle et all (I2009al ). Then, in section 
I4.2l the axial velocity and the jet shape within the ballistic portion of the jet are calculated 
through equations (|3.2j) - (|3.3p using, as initial conditions, ro = 0.5 r b , uo = B(Fr) q b /r b 
and zq(t) = z b (r) + 0.5r;,. 

4.1. Reviewing the model for r b (t) and z b (t) 
In this subsect i on we will very briefly review our model for jet formation as presented 



Gekle et al. ( 2009 ah and show its applications to predict the flow fields as well as the 



dynamics of the jet base for the impacting disc at Fr=5.1 and Fr=92 and the Worthington 
jets created after bubble pinch-off from an underwater nozzle. 

The starting point of our model is the description of the cavity collapse using a line 
of sinks on the axis of symmetry as depicted in figure 1151 After pinch-off most of this 
distribution remains intact with two notable exceptions: a hole is created between the 



bases o f the up- and downward jet and sinks accumulate around the jet base jGekle et al 
(l2009al)l. These effects are illustrated in figure (281 Based on this observation we derived 



m 



Gekle et al\ (|2009a[ ) an analytical expression for the flow potential <fi at an arbitrary 



point in the outer region (note that by construction the model is not valid inside the jet 
itself): 



2(j) 



dz' 



-oo \/r 2 + (z - z'f 



dz' 



C q b r b 



_ Zb + (z - z'f s/r 2 + (z-(z b + C Bink r b )f 



collapsing cavity 



hole 



point sink 



(4.1) 

with the order one constants C and C S i n k chosen appropriately. 

As shown in Gekle et al. ( 2009al ) this model can be used to predict the temporal 
evolution of the jet base, i.e. the widening and upwards motion of the jet base. Here 
we will restrict ourselves to show the result of this procedure for the different systems 
studied in this work, which are depicted in figure [291 

In fact, as shown in figure [30l equation (|4.ip can also be used to predict the entire 
flow field in the outer region. Figures |2"9"1I501 illustrate the rather good agreement between 
theory and numerics, which we find in all cases studied. 
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Figure 29. Time evolution of the jet base radial and axial positions, r% and zt respectively, 
taken from the simulation (magenta lines) and the analytical model (blue lines), (a) Impacting 
disc with Fr=5.1 for the upward jet (here black diamonds represent experimental data) with 
C = 4.55 and C S mk = 0.63. (b) Impacting disc with Fr=92 (upward jet) with C = 7.8 and 
Csink = 0.63. (c) Impacting disc with Fr=92 (downward jet) with C = 6.66 and C s i n k = 0.55. 
(d) Upwards jet from the underwater nozzle with C = 4.9 and C S mk = 0.76. 




Figure 30. The flow field obtained from the model with constants C = 4.55 and C s i n k = 0.63 
for the disc impacting at Fr=5.I (blue arrows) shows very good agreement with the numerically 
calculated flow field (red arrows). The region inside the jet and very close around the base is 
excluded since the model is not perfectly reliable there (due to the assumption of the observation 
point far from the base, see lGekle et al\ (|2009af )). 
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Figure 31. Comparison between the jet shape calculated using the boundary integral code 
(black line) and the one obtained integrating equations (|3.2[) - (|3.3|) using the values of ri, and z\, 
given by the model described in section |4~T1 and shown in figure [29] (red line). 



4.2. Modeling the jet breakup and drop ejection processes 

We will now take the model of the previous section one step further by combining its 
results with the analysis described in section l3.4.1l This will allow us to predict not only 
the flow field in the outer region, but also the flow inside the jet and thus the jet shape 
as a function of time. 

Once r& and Zb are obtained through the model in section 14. 1[ the axial velocity at 
the beginning of the ballistic region can be calculated as a function of known quantities 
as ito ~ B(Fr)qt,/rb, with B(Fr) the function depicted in figure [T7] Therefore, both 
the flow field and the jet shape within the ballistic region can be computed from the 
integration of equations (|3.2p - (|3.3[) using, as initial conditions, r (r) = 0.5 j-&(t), u (t) = 
B(Fi)qb(T)/rb(r) and zq(t) = z/,(t) + 0.5r(,. The comparison between the jet shape 
calculated numerically and that obtained from the model is depicted in figure [31] and 
good agreement is found. 

The capillary breakup process of the jet can be also modeled making use of our nu- 
merical results in section [3.4.21 since, through equation (|3.6|) . both Weo and sq can be 
easily expressed as a function of uq(t) = Bqb/rb, Zb(r) and ^(t), with the latter two 
functions given by the model as described above. Consequently, both the ejection and 
breakup process of the jet can be modeled with the only inputs of q c (z) and r m i n , i.e. 
quantities defined before pinch-off. 

Moreover, the trajectory of the ejected drops can also be modeled using the results of 
the previous sections. Indeed, the nearly spherical drops ejected from the tip of the jet 
follow a ballistic trajectory which can be calculated from Newton's second law as 



AttR 3 d r% 



dU d 



rop 



-p g U drop Cd 7T r eq R D 



P 5 9 



(4.2) 

3 dT 2 ry ~ arop ~ u " ' eq " u r 3 

with Udrop and Cd indicating the drop velocity and the drag coefficient, respectively. The 
drag term needs to be included since the relative variation of the drop velocity associated 
to aerodynamic effects, AUdrop/Udrop, can be estimated from equation (|4.2|) , yielding 

A Udrop Pg Udrop Cd A t flight 



Udr 



Rd 



(4.3) 



with At flight the flight time of the drop. Therefore, since Atfu g ht ~ Udrop/ 9 and Ud- 



rop 
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Uo ~ Vd B(Fr) qb/rb, the relative variation of drop velocity associated to aerodynamic 
forces is given by 

A Udr °P p a U <* r °P ^£^Fr~0(l). (4.4) 

Udrop p gR D p r\ 

It needs to be pointed out that the validity of equation (|4.2p rests on the assumption 
that drops conserve their spherical shape along their trajectory and, thus, are not atom- 
ized due to aerodynamic effects. This will be the case whenever the a erody n amic Weber 
number We a = p g Uj rop RDrdrop/o- ~ p g /pWc a < 6 [Hanson et al. (1963); Villermaux 
( 20071 120091) ]. Therefore, except for the very initial instants after cavity pinch-off, in 



which the gas Weber number could be larger than 6 - as can be inferred from figure [T4l 
-, equation (|4.2[) is valid to calculate the drop velocity. 

We are now able to calculate the trajectories of the drops ejected from the jet tip. 
Indeed, once the constants C and C s i n k of the jet base model in section rOl are properly 
chosen, they determine the values of Weo and s$j oca i, which are the only inputs for 
the model and simulations described in section I5T41 With this knowledge, the numerical 
results of the type illustrated in figures [231 and |2"T1 ( which only depend on Wes and TSq) 
allow one to calculate the size, velocity and ejection time of the first ejected droplet, the 
only inputs needed for the integration of equation 14.21 



5. Conclusions 

Using detailed boundary-integral simulations together with analytical modeling, we 
have studied the formation and breakup of the high-speed Worthington jets ejected either 
after the impact of a solid object on a liquid surface or after the pinch-off of a gas bubble 
from an underwater nozzle. To describe the phenomenon as a whole we divided the flow 
structure in two parts separated by the jet base (r^, z^): the outer region for r > r^, z < Zb 
and the jet region, extending from the jet base to the axis i.e, r < r^ and z ^ Zf,. The jet 
region is further subdivided into the three subrcgions: The axial acceleration region, where 
the radial inflow induced by the cavity collapse is decelerated radially and accelerated 
axially, the ballistic region, where fluid particles arc no longer accelerated vertically and, 
thus, conserve the axial momentum they possess at the end of the acceleration region 
and the jet tip region, which is where the jet breakup process occurs. 

We first show t hat the flow i n the outer region is well described by the analytical 
model presented in lGekle et all ( 2009al ). This model further provides a set of equations 



for the time evolution of the jet base rb(t) and Zb(t). As depicted in figures l29l and [30l 
the analytical predictions are in remarkable agreement with numerical simulations for the 
up- and downwards jets of the disc impact as well as the upwards jet created after the 
bubble pinch-off from an underwater nozzle. The model uses as its only input parameters 
the minimum radius of the cavity r min and the sink strength q c (z), both taken at the 
moment of pinch-off. 

The axial acceleration region, of characteristic length 0(rb) <C z& is where the fluid 
is decelerated in the radial direction which causes an overpressure that accelerates the 
fluid vertically. This is thus a very narrow region, localized nearby the jet base, of crucial 
importance for the jet ejection process since it is where the fluid particles transform their 
radial momentum into axial momentum. We have found the remarkable result that both 
radial (v) and axial (u) velocities, when normalized with q c (z = Zb)/rb = qb/rb, nearly 
collapse onto the same master curves for both the disc and the nozzle. Therefore, the 
values of the rcscalcd velocities (u,v)/ (qb/fb) are almost constant in time for a fixed 
value of the rescaled position r/r& < 1. We have also found that v/ (qb/rb) — for 



28 



Stephan Gekle and Jose Manuel Gordillo 



that part of jet surface whose radius is smaller than ro = 0.5 r&. Therefore, since the 
source of axial acceleration - radial deceleration of the fluid - is no longer active when 
r < ro the corresponding vertical position zq constitutes the upper boundary of the 
acceleration region. In addition, we have found that the normalized axial velocity at zo, 
u o I '(lib I V&) = B is a function which depends very weakly on time and on the Froude 
number for the impacting disc case (see figure [TT)) . 

In the slender ballistic region the axial pressure gradients arc negligible since v ~ and 
the Weber number evaluated at the beginning of the ballistic region (Wc = We^ro) is 
much larger than unity. Therefore, we have developed a ID model assuming that, in a 
first approach, fluid particles conserve their vertical velocities along the ballistic portion 
of the jet. This model allows us to calculate both the velocity field and the jet shape from 
equations (|3.2p - (|3.7p . The only input parameters are the radius, vertical position, and 
axial velocity at the beginning of the ballistic region. For the impacting disc, these values 
ro(t), zo(t), and uo(t), respectively, can be obtained directly from the analytical model of 
the outer region together with the function £?(Fr) ~ constant describing the acceleration 
region. For the underwater nozzle, the input parameters are provided directly by the 
numerical simulation. The results of this new model for the jet shape are in remarkable 
agreement with numerical simulations, as depicted in figures |2"T1 and I3T1 

Finally, we have analyzed the tip break-up region of the stretched jet. The main result 
is that the jet capillary breakup can be described as a function of two dimensionlcss 
parameters: the local Weber number Wc and the strain rate evaluated at the beginning 
of the ballistic region, ,s — du/dz(z = z ). Both quantities can again be obtained either 
from the numerical simulations or from the models of the outer and acceleration regions. 
In order to study systematically the jet breakup process as a function of these two values 
we have simulated the injection of a liquid into the atmosphere from a nozzle of constant 
radius (sec figure [4]). The real jet breakup process can then be reproduced provided that 
the values of the Weber number and the strain rate at the nozzle exit coincide with those 
at the beginning of the ballistic region, as shown in figure l23l 

We have found that the tip breakup in our physical situations is not triggered by the 
growth of perturbations coming from an external source of noise. Instead, the jet breaks 
up due to the capillary deceleration of the liquid at the tip, which produces a corrugation 
to the jet shape. Moreover, for sufficiently large values of Weo, the time evolution of 
the tip of the jet does not depend on Weo and so separately, but can be described in 
terms of the dimensionlcss parameter Weg = Wer-p Sq and the rescaled time TSq. This 
universal description allows us thus to obtain the size of the droplet ejected from the 
tip (cf. figure I2T[) if Weo and sq are known fr om either simulatio ns, measurements, or 
analytical models such as the one described in Gekle et all ( 2009a ). 

In summary, our description of Worthington jets created by the impact of a solid object 
on a liquid surface allows us to predict the jet base dynamics, the jet shape, and even 
the ejection of drops from the tip of the jet based only on the knowledge of the minimum 
radius of the cavity before the jet emerges and the sink distribution at pinch-off. 
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